## estimates an individual-level regession of 2016 voting against 2012 voting plus controls

library(foreign)
library(lmtest)
library(effects)
library(nnet)
library(MASS)
library(tidyverse)
library(plyr) 
library(dplyr)
library(haven)
library(simcf)
library(ggplot2)
library(multiwayvcov)
library(stargazer)
options(scipen=20)



## attach the vbf object unless it is already attached
if (!length(grep("vbf_output", search()))) {
  cat ("Attaching vbf... this may take a few moments...")
  attach("vbf_output")
}

## predicting 2016 turnout ##

# model 

# get the factor vars ready with appropriate base categories
vbf$race <- as.factor(vbf$race)
vbf$age_category <- as.factor(vbf$age_category)
vbf <- within(vbf, race <- relevel(race, ref = "white"), age_category <- relevel(age_category, ref = "age_middle"))

# Function to calculate clustered SEs 
robustse.f <- function(model, cluster, df_correction) {
  require(sandwich)
  require(lmtest)
  require(multiwayvcov)
  if(missing(cluster)) {
    name <- deparse(substitute(model))
    modelname <- paste(name,"rob",sep=".")
    model$se <- coeftest(model, vcov=vcovHC(model,"HC1"))[,2]
    model$vcovHC <- vcovHC(model,"HC1")
    assign(modelname,model,envir = .GlobalEnv)
    coeftest(model, vcov=vcovHC(model,"HC1"))
  } else {
    name <- deparse(substitute(model))
    modelname <- paste(name,"rob",sep=".")
    vcovCL <- cluster.vcov(model, cluster, df_correction = df_correction)
    model$vcovCL <- vcovCL
    modelname <- paste(name,"clustrob",sep=".")
    model$se <- coeftest(model, vcovCL)[,2]
    assign(modelname,model,envir = .GlobalEnv)
    coeftest(model, vcovCL)
  }
}

# all institutional changes in one model 

summary(model_allchanges <- glm(factor(outcome2016 != "novote2016") ~ totalhoursChg + evehoursChg + 
                                     sathoursChg+sunhoursChg+sitesChg+
                                     factor(race) + factor(female) + factor(party) + 
                                     factor(age_category) + factor(outcome2012 != "novote2012"),
                                   data = vbf %>% filter (voter_status_desc_2012 == "ACTIVE",
                                                          voter_status_desc_2016 == "ACTIVE",
                                                          ethnicity == "latino_not_yes",
                                                          (race == "black" | race == "white"),
                                                          (party == "dem" | party=="rep")),
                                   family = "binomial"))


model_totalhoursChg_cluster <- robustse.f(model_allchanges, model_allchanges$data$COUNTY_NAM, df_correction = TRUE)


stargazer(model_totalhoursChg_cluster,
          no.space = TRUE, 
          title = "Changes in EIP availability and 2016 voter turnout",
          dep.var.labels.include = TRUE,
          covariate.labels=c("Percent Change: Total Hours",
                             "Percent Change: Evening Hours",
                             "Percent Change: Saturday Hours",
                             "Percent Change: Sunday Hours",
                             "Percent Change: EIP Sites",
                             "Black","Female","Republican",
                             "Age: 18-34","Age: 65+","Voted in 2012",
                             "Constant"),
          order=c(1,2,3,4,5,6,7,8,10,9,11,12),
          label = "tab:base_regression",
          model.numbers = FALSE,
          star.cutoffs = c(0.05, 0.01, 0.001),
          out = "models_Chg_cluster_base.tex",
          dep.var.labels = "")

# interactions #

### THREEWAY INTERACTIONS FOR OTHER INSTITUTIONAL CHANGES

# total hours
summary(model_totalhoursChg <- glm(factor(outcome2016 != "novote2016") ~ totalhoursChg * factor(race) * factor(party) + 
                                     evehoursChg + sathoursChg+sunhoursChg+sitesChg+ factor(female) + 
                                     factor(age_category) + factor(outcome2012 != "novote2012"),
                                   data = vbf %>% filter (voter_status_desc_2012 == "ACTIVE",
                                                          voter_status_desc_2016 == "ACTIVE",
                                                          ethnicity == "latino_not_yes",
                                                          (race == "black" | race == "white"),
                                                          (party == "dem" | party=="rep")),
                                   family = "binomial"))


model_totalhoursChg_interact_cluster <- robustse.f(model_totalhoursChg, model_totalhoursChg$data$COUNTY_NAM, df_correction = TRUE)

stargazer(model_totalhoursChg_interact_cluster,
          no.space = TRUE, 
          title = "Changes in total EIP hours, race, party and 2016 voter turnout",
          dep.var.labels.include = TRUE,
          order=c(12,13,14,15,1,2,3,4,5,6,7,8,9,10,11),
          covariate.labels=c("Total Hours*Black","Total Hours*Republican",
                             "Black*Republican","Total Hours*Black*Republican",
                             "Total Hours","Black","Republican",
                             "Evening Hours","Saturday Hours","Sunday Hours","Sites",
                             "Female","Age: 18-34","Age: 65+","Voted 2012","Constant"),
          label = "tab:totalhours_interact",
          model.numbers = FALSE,
          star.cutoffs = c(0.05, 0.01, 0.001),
          out = "totalhoursChg_cluster_interact.tex",
          dep.var.labels = "")

# evening hours 
summary(model_evehoursChg <- glm(factor(outcome2016 != "novote2016") ~ evehoursChg * factor(race) * factor(party) + 
                                   totalhoursChg+sathoursChg+sunhoursChg+sitesChg+
                                   factor(female) + factor(age_category) + factor(outcome2012 != "novote2012"),
                                 data = vbf %>% filter (voter_status_desc_2012 == "ACTIVE",
                                                        voter_status_desc_2016 == "ACTIVE",
                                                        ethnicity == "latino_not_yes",
                                                        (race == "black" | race == "white"),
                                                        (party == "dem" | party=="rep")),
                                 family = "binomial"))

model_evehoursChg_interact_cluster <- robustse.f(model_evehoursChg, model_evehoursChg$data$COUNTY_NAM, df_correction = TRUE)

stargazer(model_evehoursChg_interact_cluster,
          no.space = TRUE, 
          title = "Changes in evening EIP hours, race, party and 2016 voter turnout",
          dep.var.labels.include = TRUE,
          order=c(12,13,14,15,1,2,3,4,5,6,7,8,9,10,11),
          covariate.labels=c("Evening Hours*Black","Evening Hours*Republican",
                             "Black*Republican","Evening Hours*Black*Republican",
                             "Evening Hours","Black","Republican",
                             "Total Hours","Saturday Hours","Sunday Hours","Sites",
                             "Female","Age: 18-34","Age: 65+","Voted 2012","Constant"),
          label = "tab:evehours_interact",
          model.numbers = FALSE,
          star.cutoffs = c(0.05, 0.01, 0.001),
          out = "evehoursChg_cluster_interact.tex",
          dep.var.labels = "")



# saturday hours

summary(model_sathoursChg <- glm(factor(outcome2016 != "novote2016") ~ sathoursChg * factor(race) * factor(party) +
                                   totalhoursChg+evehoursChg+sunhoursChg+sitesChg+factor(female) + 
                                   factor(age_category) + factor(outcome2012 != "novote2012"),
                                 data = vbf %>% filter (voter_status_desc_2012 == "ACTIVE",
                                                        voter_status_desc_2016 == "ACTIVE",
                                                        ethnicity == "latino_not_yes",
                                                        (race == "black" | race == "white"),
                                                        (party == "dem" | party=="rep")),
                                 family = "binomial"))

model_sathoursChg_interact_cluster <- robustse.f(model_sathoursChg, model_sathoursChg$data$COUNTY_NAM, df_correction = TRUE)

stargazer(model_sathoursChg_interact_cluster,
          no.space = TRUE, 
          title = "Changes in Saturday EIP hours, race, party and 2016 voter turnout",
          dep.var.labels.include = TRUE,
          order=c(12,13,14,15,1,2,3,4,5,6,7,8,9,10,11),
          covariate.labels=c("Saturday Hours*Black","Saturday Hours*Republican",
                             "Black*Republican","Saturday Hours*Black*Republican",
                             "Saturday Hours","Black","Republican",
                             "Total Hours","Evening Hours","Sunday Hours","Sites",
                             "Female", "Age: 18-34","Age: 65+","Voted 2012","Constant"),
          label = "tab:satuhours_interact",
          model.numbers = FALSE,
          star.cutoffs = c(0.05, 0.01, 0.001),
          out = "sathoursChg_cluster_interact.tex",
          dep.var.labels = "")


# sunday hours

summary(model_sunhoursChg <- glm(factor(outcome2016 != "novote2016") ~ sunhoursChg * factor(race) * factor(party) + 
                                   totalhoursChg+evehoursChg+sathoursChg+sitesChg+factor(female) + 
                                   factor(age_category) + factor(outcome2012 != "novote2012"),
                                 data = vbf %>% filter (voter_status_desc_2012 == "ACTIVE",
                                                        voter_status_desc_2016 == "ACTIVE",
                                                        ethnicity == "latino_not_yes",
                                                        (race == "black" | race == "white"),
                                                        (party == "dem" | party=="rep")),
                                 family = "binomial"))

model_sunhoursChg_interact_cluster <- robustse.f(model_sunhoursChg, model_sunhoursChg$data$COUNTY_NAM, df_correction = TRUE)

stargazer(model_sunhoursChg_interact_cluster,
          no.space = TRUE, 
          title = "Changes in Sunday EIP hours, race, party and 2016 voter turnout",
          dep.var.labels.include = TRUE,
          order=c(12,13,14,15,1,2,3,4,5,6,7,8,9,10,11),
          covariate.labels=c("Sunday Hours*Black","Sunday Hours*Republican",
                             "Black*Republican","Sunday Hours*Black*Republican",
                             "Sunday Hours","Black","Republican",
                             "Total Hours","Evening Hours","Saturday Hours","Sites",
                             "Female","Age: 18-34","Age: 65+","Voted 2012","Constant"),
          label = "tab:sunhours_interact",
          model.numbers = FALSE,
          star.cutoffs = c(0.05, 0.01, 0.001),
          out = "sunhoursChg_cluster_interact.tex",
          dep.var.labels = "")

# sites

summary(model_sitesChg <- glm(factor(outcome2016 != "novote2016") ~ sitesChg * factor(race) * factor(party) + 
                                totalhoursChg+evehoursChg+sathoursChg+sunhoursChg+
                                factor(female) + factor(age_category) + factor(outcome2012 != "novote2012"),
                              data = vbf %>% filter (voter_status_desc_2012 == "ACTIVE",
                                                     voter_status_desc_2016 == "ACTIVE",
                                                     ethnicity == "latino_not_yes",
                                                     (race == "black" | race == "white"),
                                                     (party == "dem" | party=="rep")),
                              family = "binomial"))

model_sitesChg_interact_cluster <- robustse.f(model_sitesChg, model_sitesChg$data$COUNTY_NAM, df_correction = TRUE)

stargazer(model_sitesChg_interact_cluster,
          no.space = TRUE, 
          title = "Changes in EIP sites, race, party and 2016 voter turnout",
          dep.var.labels.include = TRUE,
          order=c(12,13,14,15,1,2,3,4,5,6,7,8,9,10,11),
          covariate.labels=c("Sites Hours*Black","Sites Hours*Republican",
                             "Black*Republican","Sites Hours*Black*Republican",
                             "Sites Hours","Black","Republican",
                             "Total Hours","Evening Hours","Saturday Hours","Sunday Hours",
                             "Female", "Age: 18-34","Age: 65+","Voted 2012","Constant"),
          label = "tab:sites_interact",
          model.numbers = FALSE,
          star.cutoffs = c(0.05, 0.01, 0.001),
          out = "sitesChg_cluster_interact.tex",
          dep.var.labels = "")

stop()



